1 | #pragma rtGlobals=1 // Use modern global access method. |
---|
2 | #pragma IgorVersion = 6.0 |
---|
3 | |
---|
4 | // currently, there is NO XOP version of RPA, since there is an extra input parameter |
---|
5 | // wave that must be carried into the function. |
---|
6 | // this can be done with the STRUCT, if extra space is allocated for such |
---|
7 | // |
---|
8 | // |
---|
9 | Proc PlotRPAForm(num,qmin,qmax,nCASE) |
---|
10 | Variable num=100,qmin=0.001,qmax=0.5 |
---|
11 | Variable/g gCASE //Global Variable |
---|
12 | Variable nCASE |
---|
13 | Prompt num "Enter number of data points: " |
---|
14 | Prompt qmin "Enter minimum Q value (A^-1): " |
---|
15 | Prompt qmax "Enter maximum Q value (A^-1): " |
---|
16 | Prompt nCASE, "Choose one of the following cases:",popup "CASE 1: C/D BINARY MIXTURE OF HOMOPOLYMERS;" |
---|
17 | "CASE 2: C-D DIBLOCK COPOLYMER;" |
---|
18 | "CASE 3: B/C/D TERNARY MIXTURE OF HOMOPOLYMERS;" |
---|
19 | "CASE 4: B/C-D MIXTURE OF HOMOPOLYMER B AND DIBLOCK COPOLYMER C-D;" |
---|
20 | "CASE 5: B-C-D TRIBLOCK COPOLYMER;" |
---|
21 | "CASE 6: A/B/C/D QUATERNARY MIXTURE OF HOMOPOLYMERS;" |
---|
22 | "CASE 7: A/B/C-D MIXTURE OF TWO HOMOPOLYMERS A/B AND A DIBLOCK C-D;" |
---|
23 | "CASE 8: A/B-C-D MIXTURE OF A HOMOPOLYMER A AND A TRIBLOCK B-C-D;" |
---|
24 | "CASE 9: A-B/C-D MIXTURE OF TWO DIBLOCK COPOLYMERS A-B AND C-D;" |
---|
25 | "CASE 10: A-B-C-D FOUR-BLOCK COPOLYMER" |
---|
26 | |
---|
27 | Make/O/D/n=(num) xwave_rpa,ywave_rpa |
---|
28 | //xwave_rpa=qmin+x*(qmax-qmin)/num |
---|
29 | //switch to log-scaling of the q-values |
---|
30 | xwave_rpa=alog(log(qmin) + x*((log(qmax)-log(qmin))/num)) |
---|
31 | |
---|
32 | gCASE = nCASE |
---|
33 | // print gCASE |
---|
34 | |
---|
35 | IF(gCASE <=2) |
---|
36 | Make/O/D inputvalues={1000,0.5,100,1.e-12,1000,0.5,100,0.e-12} |
---|
37 | make/o/t inputnames={"Deg. Polym. Nc","Vol. Frac. Phic","Spec. Vol. Vc","Scatt. Length Lc","Deg. Polym. Nd","Vol. Frac. Phid","Spec. Vol. Vd","Scatt. LengthLd"} |
---|
38 | Edit inputnames,inputvalues |
---|
39 | Make/O/D coef_rpa = {5,5,-.0004,1,0} |
---|
40 | make/o/t parameters_rpa = {"Seg. Length bc","Seg. Length bd","Chi Param. Kcd","scale","Background"} |
---|
41 | Edit parameters_rpa,coef_rpa |
---|
42 | ENDIF |
---|
43 | |
---|
44 | IF((gCASE >2) %& (gCASE<=5)) |
---|
45 | Make/O/D inputvalues={1000,0.33,100,1.e-12,1000,0.33,100,1.e-12,1000,0.33,100,0.e-12} |
---|
46 | make/o/t inputnames={"Deg. Polym. Nb","Vol. Frac. Phib","Spec. Vol. Vb","Scatt. Length Lb","Deg. Polym. Nc","Vol. Frac. Phic","Spec. Vol. Vc","Scatt. Length Lc","Deg. Polym. Nd","Vol. Frac. Phid","Spec. Vol. Vd","Scatt. Length Ld"} |
---|
47 | Edit inputnames,inputvalues |
---|
48 | Make/O/D coef_rpa = {5,5,5,-.0004,-.0004,-.0004,1,0} |
---|
49 | make/o/t parameters_rpa = {"Seg. Length bb","Seg. Length bc","Seg. Length bd","Chi Param. Kbc","Chi Param. Kbd","Chi Param. Kcd","scale","Background"} |
---|
50 | Edit parameters_rpa,coef_rpa |
---|
51 | ENDIF |
---|
52 | |
---|
53 | IF(gCASE >5) |
---|
54 | Make/O/D inputvalues={1000,0.25,100,1.e-12,1000,0.25,100,1.e-12,1000,0.25,100,1.e-12,1000,0.25,100,0.e-12} |
---|
55 | make/o/t inputnames={"Deg. Polym. Na","Vol. Frac. Phia","Spec. Vol. Va","Scatt. Length La","Deg. Polym. Nb","Vol. Frac. Phib","Spec. Vol. Vb","Scatt. Length Lb","Deg. Polym. Nc","Vol. Frac. Phic","Spec. Vol. Vc","Scatt. Length Lc","Deg. Polym. Nd","Vol. Frac. Phid","Spec. Vol. Vd","Scatt. Length Ld"} |
---|
56 | Edit inputnames,inputvalues |
---|
57 | Make/O/D coef_rpa = {5,5,5,5,-.0004,-.0004,-.0004,-.0004,-.0004,-.0004,1,0} |
---|
58 | make/o/t parameters_rpa = {"Seg. Length ba","Seg. Length bb","Seg. Length bc","Seg. Length bd","Chi Param. Kab","Chi Param. Kac","Chi Param. Kad","Chi Param. Kbc","Chi Param. Kbd","Chi Param. Kcd","scale","Background"} |
---|
59 | Edit parameters_rpa,coef_rpa |
---|
60 | ENDIF |
---|
61 | |
---|
62 | Variable/G root:g_rpa |
---|
63 | g_rpa := RPAForm(coef_rpa,ywave_rpa,xwave_rpa) |
---|
64 | Display ywave_rpa vs xwave_rpa |
---|
65 | // ModifyGraph log=1,marker=29,msize=2,mode=4 |
---|
66 | Label bottom "Q (A\\S-1\\M)" |
---|
67 | Label left "Intensity (cm\\S-1\\M)" |
---|
68 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
69 | End |
---|
70 | |
---|
71 | // - sets up a dependency to a wrapper, not the actual SmearedModelFunction |
---|
72 | Proc PlotSmearedRPAForm(str,nCASE) |
---|
73 | String str |
---|
74 | Variable nCASE |
---|
75 | Prompt str,"Pick the data folder containing the resolution you want",popup,getAList(4) |
---|
76 | Prompt nCASE, "Choose one of the following cases:",popup "CASE 1: C/D BINARY MIXTURE OF HOMOPOLYMERS;" |
---|
77 | "CASE 2: C-D DIBLOCK COPOLYMER;" |
---|
78 | "CASE 3: B/C/D TERNARY MIXTURE OF HOMOPOLYMERS;" |
---|
79 | "CASE 4: B/C-D MIXTURE OF HOMOPOLYMER B AND DIBLOCK COPOLYMER C-D;" |
---|
80 | "CASE 5: B-C-D TRIBLOCK COPOLYMER;" |
---|
81 | "CASE 6: A/B/C/D QUATERNARY MIXTURE OF HOMOPOLYMERS;" |
---|
82 | "CASE 7: A/B/C-D MIXTURE OF TWO HOMOPOLYMERS A/B AND A DIBLOCK C-D;" |
---|
83 | "CASE 8: A/B-C-D MIXTURE OF A HOMOPOLYMER A AND A TRIBLOCK B-C-D;" |
---|
84 | "CASE 9: A-B/C-D MIXTURE OF TWO DIBLOCK COPOLYMERS A-B AND C-D;" |
---|
85 | "CASE 10: A-B-C-D FOUR-BLOCK COPOLYMER" |
---|
86 | |
---|
87 | SetDataFolder $("root:"+str) |
---|
88 | |
---|
89 | Variable/g gCASE //Global Variable |
---|
90 | |
---|
91 | // if any of the resolution waves are missing => abort |
---|
92 | if(ResolutionWavesMissingDF(str)) //updated to NOT use global strings (in GaussUtils) |
---|
93 | Abort |
---|
94 | endif |
---|
95 | |
---|
96 | gCASE = nCASE |
---|
97 | // print gCASE |
---|
98 | |
---|
99 | |
---|
100 | IF(gCASE <=2) |
---|
101 | Make/O/D inputvalues={1000,0.5,100,1.e-12,1000,0.5,100,0.e-12} |
---|
102 | make/o/t inputnames={"Deg. Polym. Nc","Vol. Frac. Phic","Spec. Vol. Vc","Scatt. Length Lc","Deg. Polym. Nd","Vol. Frac. Phid","Spec. Vol. Vd","Scatt. LengthLd"} |
---|
103 | Edit inputnames,inputvalues |
---|
104 | Make/O/D smear_coef_rpa = {5,5,-.0004,1,0} |
---|
105 | make/o/t smear_parameters_rpa = {"Seg. Length bc","Seg. Length bd","Chi Param. Kcd","scale","Background"} |
---|
106 | Edit smear_parameters_rpa,smear_coef_rpa |
---|
107 | ENDIF |
---|
108 | |
---|
109 | IF((gCASE >2) %& (gCASE<=5)) |
---|
110 | Make/O/D inputvalues={1000,0.33,100,1.e-12,1000,0.33,100,1.e-12,1000,0.33,100,0.e-12} |
---|
111 | make/o/t inputnames={"Deg. Polym. Nb","Vol. Frac. Phib","Spec. Vol. Vb","Scatt. Length Lb","Deg. Polym. Nc","Vol. Frac. Phic","Spec. Vol. Vc","Scatt. Length Lc","Deg. Polym. Nd","Vol. Frac. Phid","Spec. Vol. Vd","Scatt. Length Ld"} |
---|
112 | Edit inputnames,inputvalues |
---|
113 | Make/O/D smear_coef_rpa = {5,5,5,-.0004,-.0004,-.0004,1,0} |
---|
114 | make/o/t smear_parameters_rpa = {"Seg. Length bb","Seg. Length bc","Seg. Length bd","Chi Param. Kbc","Chi Param. Kbd","Chi Param. Kcd","scale","Background"} |
---|
115 | Edit smear_parameters_rpa,smear_coef_rpa |
---|
116 | ENDIF |
---|
117 | |
---|
118 | IF(gCASE >5) |
---|
119 | Make/O/D inputvalues={1000,0.25,100,1.e-12,1000,0.25,100,1.e-12,1000,0.25,100,1.e-12,1000,0.25,100,0.e-12} |
---|
120 | make/o/t inputnames={"Deg. Polym. Na","Vol. Frac. Phia","Spec. Vol. Va","Scatt. Length La","Deg. Polym. Nb","Vol. Frac. Phib","Spec. Vol. Vb","Scatt. Length Lb","Deg. Polym. Nc","Vol. Frac. Phic","Spec. Vol. Vc","Scatt. Length Lc","Deg. Polym. Nd","Vol. Frac. Phid","Spec. Vol. Vd","Scatt. Length Ld"} |
---|
121 | Edit inputnames,inputvalues |
---|
122 | Make/O/D smear_coef_rpa = {5,5,5,5,-.0004,-.0004,-.0004,-.0004,-.0004,-.0004,1,0} |
---|
123 | make/o/t smear_parameters_rpa = {"Seg. Length ba","Seg. Length bb","Seg. Length bc","Seg. Length bd","Chi Param. Kab","Chi Param. Kac","Chi Param. Kad","Chi Param. Kbc","Chi Param. Kbd","Chi Param. Kcd","scale","Background"} |
---|
124 | Edit smear_parameters_rpa,smear_coef_rpa |
---|
125 | ENDIF |
---|
126 | |
---|
127 | // output smeared intensity wave, dimensions are identical to experimental QSIG values |
---|
128 | // make extra copy of experimental q-values for easy plotting |
---|
129 | Duplicate/O $(str+"_q") smeared_rpa,smeared_qvals |
---|
130 | SetScale d,0,0,"1/cm",smeared_rpa |
---|
131 | |
---|
132 | Variable/G gs_rpa=0 |
---|
133 | gs_rpa := fSmearedRPAForm(smear_coef_rpa,smeared_rpa,smeared_qvals) //this wrapper fills the STRUCT |
---|
134 | |
---|
135 | Display smeared_rpa vs smeared_qvals |
---|
136 | ModifyGraph log=1,marker=29,msize=2,mode=4 |
---|
137 | Label bottom "q (A\\S-1\\M)" |
---|
138 | Label left "Intensity (cm\\S-1\\M)" |
---|
139 | AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) |
---|
140 | |
---|
141 | SetDataFolder root: |
---|
142 | End |
---|
143 | |
---|
144 | /////////////////////////////////////////////////////////////// |
---|
145 | //AAO version, uses XOP if available |
---|
146 | // simply calls the original single point calculation with |
---|
147 | // a wave assignment (this will behave nicely if given point ranges) |
---|
148 | Function RPAForm(cw,yw,xw) : FitFunc |
---|
149 | Wave cw,yw,xw |
---|
150 | |
---|
151 | #if exists("RPAFormX") |
---|
152 | yw = RPAFormX(cw,xw) |
---|
153 | #else |
---|
154 | yw = fRPAForm(cw,xw) |
---|
155 | #endif |
---|
156 | return(0) |
---|
157 | End |
---|
158 | |
---|
159 | |
---|
160 | Function fRPAForm(w,x) : FitFunc |
---|
161 | Wave w |
---|
162 | Variable x |
---|
163 | |
---|
164 | Wave var=$"inputvalues" |
---|
165 | Nvar lCASE=gCASE |
---|
166 | // print lCASE |
---|
167 | // Variable lCASE |
---|
168 | // RANDOM PHASE APPROXIMATION FOR A FOUR-BLOCK COPOLYMER A-B-C-D |
---|
169 | // THIS FORMALISM APPLIES TO MULTICOMPONENT POLYMER MIXTURES |
---|
170 | // IN THE HOMOGENEOUS (MIXED) PHASE REGION ONLY. |
---|
171 | |
---|
172 | // THIS GENERAL CASE INCLUDES A MULTITUDE OF (TEN) SPECIAL CASES. |
---|
173 | |
---|
174 | // HERE ARE THE VARIOUS CASES COVERED: |
---|
175 | // CASE 1: C/D BINARY MIXTURE OF HOMOPOLYMERS |
---|
176 | // CASE 2: C-D DIBLOCK COPOLYMER |
---|
177 | // CASE 3: B/C/D TERNARY MIXTURE OF HOMOPOLYMERS |
---|
178 | // CASE 4: B/C-D MIXTURE OF HOMOPOLYMER B AND DIBLOCK COPOLYMER C-D |
---|
179 | // CASE 5: B-C-D TRIBLOCK COPOLYMER |
---|
180 | // CASE 6: A/B/C/D QUATERNARY MIXTURE OF HOMOPOLYMERS |
---|
181 | // CASE 7: A/B/C-D MIXTURE OF TWO HOMOPOLYMERS A/B AND A DIBLOCK C-D |
---|
182 | // CASE 8: A/B-C-D MIXTURE OF A HOMOPOLYMER A AND A TRIBLOCK B-C-D |
---|
183 | // CASE 9: A-B/C-D MIXTURE OF TWO DIBLOCK COPOLYMERS A-B AND C-D |
---|
184 | // CASE 10: A-B-C-D FOUR-BLOCK COPOLYMER |
---|
185 | |
---|
186 | // B. HAMMOUDA, NIST, JULY 1998 |
---|
187 | |
---|
188 | Variable Na,Nb,Nc,Nd,Nab,Nac,Nad,Nba,Nbc,Nbd,Nca,Ncb,Ncd |
---|
189 | Variable Phia,Phib,Phic,Phid,Phiab,Phiac,Phiad |
---|
190 | Variable Phiba,Phibc,Phibd,Phica,Phicb,Phicd,Phida,Phidb,Phidc |
---|
191 | Variable va,vb,vc,vd,vab,vac,vad,vba,vbc,vbd,vca,vcb,vcd,vda,vdb,vdc |
---|
192 | Variable m |
---|
193 | Variable ba,bb,bc,bd |
---|
194 | Variable Q |
---|
195 | Variable Xa,Xb,Xc,Xd |
---|
196 | Variable Paa,S0aa,Pab,S0ab,Pac,S0ac,Pad,S0ad |
---|
197 | Variable Pba,S0ba,Pbb,S0bb,Pbc,S0bc,Pbd,S0bd |
---|
198 | Variable Pca,S0ca,Pcb,S0cb,Pcc,S0cc,Pcd,S0cd |
---|
199 | Variable Pda,S0da,Pdb,S0db,Pdc,S0dc,Pdd,S0dd |
---|
200 | Variable Kaa,Kab,Kac,Kad,Kba,Kbb,Kbc,Kbd |
---|
201 | Variable Kca,Kcb,Kcc,Kcd,Kda,Kdb,Kdc,Kdd |
---|
202 | Variable Zaa,Zab,Zac,Zba,Zbb,Zbc,Zca,Zcb,Zcc |
---|
203 | Variable DenT,T11,T12,T13,T21,T22,T23,T31,T32,T33 |
---|
204 | Variable Y1,Y2,Y3,X11,X12,X13,X21,X22,X23,X31,X32,X33 |
---|
205 | Variable ZZ,DenQ1,DenQ2,DenQ3,DenQ,Q11,Q12,Q13,Q21,Q22,Q23,Q31,Q32,Q33 |
---|
206 | Variable N11,N12,N13,N21,N22,N23,N31,N32,N33 |
---|
207 | Variable M11,M12,M13,M21,M22,M23,M31,M32,M33 |
---|
208 | Variable S11,S12,S13,S14,S21,S22,S23,S24 |
---|
209 | Variable S31,S32,S33,S34,S41,S42,S43,S44 |
---|
210 | Variable La,Lb,Lc,Ld,Lad,Lbd,Lcd,Nav,Int |
---|
211 | Variable scale |
---|
212 | Variable Background |
---|
213 | Variable ii=0 //to fool curve fitting dialog to let you pick the coefficient wave |
---|
214 | |
---|
215 | Na=1000 |
---|
216 | Nb=1000 |
---|
217 | Nc=1000 |
---|
218 | Nd=1000 //DEGREE OF POLYMERIZATION |
---|
219 | Phia=0.25 |
---|
220 | Phib=0.25 |
---|
221 | Phic=0.25 |
---|
222 | Phid=0.25 //VOL FRACTION |
---|
223 | Kab=-.0004 |
---|
224 | Kac=-.0004 |
---|
225 | Kad=-.0004 //CHI PARAM |
---|
226 | Kbc=-.0004 |
---|
227 | Kbd=-.0004 |
---|
228 | Kcd=-.0004 |
---|
229 | La=1.E-12 |
---|
230 | Lb=1.E-12 |
---|
231 | Lc=1.E-12 |
---|
232 | Ld=0.E-12 //SCATT. LENGTH |
---|
233 | va=100 |
---|
234 | vb=100 |
---|
235 | vc=100 |
---|
236 | vd=100 //SPECIFIC VOLUME |
---|
237 | ba=5 |
---|
238 | bb=5 |
---|
239 | bc=5 |
---|
240 | bd=5 //SEGMENT LENGTH |
---|
241 | |
---|
242 | IF (lCASE <= 2) |
---|
243 | Phia=0.0000001 |
---|
244 | Phib=0.0000001 |
---|
245 | Phic=0.5 |
---|
246 | Phid=0.5 |
---|
247 | Nc=var[0] |
---|
248 | Phic=var[1] |
---|
249 | vc=var[2] |
---|
250 | Lc=var[3] |
---|
251 | bc=w[ii+0] |
---|
252 | Nd=var[4] |
---|
253 | Phid=var[5] |
---|
254 | vd=var[6] |
---|
255 | Ld=var[7] |
---|
256 | bd=w[ii+1] |
---|
257 | Kcd=w[ii+2] |
---|
258 | scale=w[ii+3] |
---|
259 | background=w[ii+4] |
---|
260 | ENDIF |
---|
261 | |
---|
262 | IF ((lCASE > 2) %& (lCASE <= 5)) |
---|
263 | Phia=0.0000001 |
---|
264 | Phib=0.333333 |
---|
265 | Phic=0.333333 |
---|
266 | Phid=0.333333 |
---|
267 | Nb=var[0] |
---|
268 | Phib=var[1] |
---|
269 | vb=var[2] |
---|
270 | Lb=var[3] |
---|
271 | bb=w[ii+0] |
---|
272 | Nc=var[4] |
---|
273 | Phic=var[5] |
---|
274 | vc=var[6] |
---|
275 | Lc=var[7] |
---|
276 | bc=w[ii+1] |
---|
277 | Nd=var[8] |
---|
278 | Phid=var[9] |
---|
279 | vd=var[10] |
---|
280 | Ld=var[11] |
---|
281 | bd=w[ii+2] |
---|
282 | Kbc=w[ii+3] |
---|
283 | Kbd=w[ii+4] |
---|
284 | Kcd=w[ii+5] |
---|
285 | scale=w[ii+6] |
---|
286 | background=w[ii+7] |
---|
287 | ENDIF |
---|
288 | |
---|
289 | IF (lCASE > 5) |
---|
290 | Phia=0.25 |
---|
291 | Phib=0.25 |
---|
292 | Phic=0.25 |
---|
293 | Phid=0.25 |
---|
294 | Na=var[0] |
---|
295 | Phia=var[1] |
---|
296 | va=var[2] |
---|
297 | La=var[3] |
---|
298 | ba=w[ii+0] |
---|
299 | Nb=var[4] |
---|
300 | Phib=var[5] |
---|
301 | vb=var[6] |
---|
302 | Lb=var[7] |
---|
303 | bb=w[ii+1] |
---|
304 | Nc=var[8] |
---|
305 | Phic=var[9] |
---|
306 | vc=var[10] |
---|
307 | Lc=var[11] |
---|
308 | bc=w[ii+2] |
---|
309 | Nd=var[12] |
---|
310 | Phid=var[13] |
---|
311 | vd=var[14] |
---|
312 | Ld=var[15] |
---|
313 | bd=w[ii+3] |
---|
314 | Kab=w[ii+4] |
---|
315 | Kac=w[ii+5] |
---|
316 | Kad=w[ii+6] |
---|
317 | Kbc=w[ii+7] |
---|
318 | Kbd=w[ii+8] |
---|
319 | Kcd=w[ii+9] |
---|
320 | scale=w[ii+10] |
---|
321 | background=w[ii+11] |
---|
322 | ENDIF |
---|
323 | |
---|
324 | Nab=(Na*Nb)^(0.5) |
---|
325 | Nac=(Na*Nc)^(0.5) |
---|
326 | Nad=(Na*Nd)^(0.5) |
---|
327 | Nbc=(Nb*Nc)^(0.5) |
---|
328 | Nbd=(Nb*Nd)^(0.5) |
---|
329 | Ncd=(Nc*Nd)^(0.5) |
---|
330 | |
---|
331 | vab=(va*vb)^(0.5) |
---|
332 | vac=(va*vc)^(0.5) |
---|
333 | vad=(va*vd)^(0.5) |
---|
334 | vbc=(vb*vc)^(0.5) |
---|
335 | vbd=(vb*vd)^(0.5) |
---|
336 | vcd=(vc*vd)^(0.5) |
---|
337 | |
---|
338 | Phiab=(Phia*Phib)^(0.5) |
---|
339 | Phiac=(Phia*Phic)^(0.5) |
---|
340 | Phiad=(Phia*Phid)^(0.5) |
---|
341 | Phibc=(Phib*Phic)^(0.5) |
---|
342 | Phibd=(Phib*Phid)^(0.5) |
---|
343 | Phicd=(Phic*Phid)^(0.5) |
---|
344 | |
---|
345 | Q=x |
---|
346 | Xa=Q^2*ba^2*Na/6 |
---|
347 | Xb=Q^2*bb^2*Nb/6 |
---|
348 | Xc=Q^2*bc^2*Nc/6 |
---|
349 | Xd=Q^2*bd^2*Nd/6 |
---|
350 | |
---|
351 | Paa=2*(Exp(-Xa)-1+Xa)/Xa^2 |
---|
352 | S0aa=Na*Phia*va*Paa |
---|
353 | Pab=((1-Exp(-Xa))/Xa)*((1-Exp(-Xb))/Xb) |
---|
354 | S0ab=(Phiab*vab*Nab)*Pab |
---|
355 | Pac=((1-Exp(-Xa))/Xa)*Exp(-Xb)*((1-Exp(-Xc))/Xc) |
---|
356 | S0ac=(Phiac*vac*Nac)*Pac |
---|
357 | Pad=((1-Exp(-Xa))/Xa)*Exp(-Xb-Xc)*((1-Exp(-Xd))/Xd) |
---|
358 | S0ad=(Phiad*vad*Nad)*Pad |
---|
359 | |
---|
360 | S0ba=S0ab |
---|
361 | Pbb=2*(Exp(-Xb)-1+Xb)/Xb^2 |
---|
362 | S0bb=Nb*Phib*vb*Pbb |
---|
363 | Pbc=((1-Exp(-Xb))/Xb)*((1-Exp(-Xc))/Xc) |
---|
364 | S0bc=(Phibc*vbc*Nbc)*Pbc |
---|
365 | Pbd=((1-Exp(-Xb))/Xb)*Exp(-Xc)*((1-Exp(-Xd))/Xd) |
---|
366 | S0bd=(Phibd*vbd*Nbd)*Pbd |
---|
367 | |
---|
368 | S0ca=S0ac |
---|
369 | S0cb=S0bc |
---|
370 | Pcc=2*(Exp(-Xc)-1+Xc)/Xc^2 |
---|
371 | S0cc=Nc*Phic*vc*Pcc |
---|
372 | Pcd=((1-Exp(-Xc))/Xc)*((1-Exp(-Xd))/Xd) |
---|
373 | S0cd=(Phicd*vcd*Ncd)*Pcd |
---|
374 | |
---|
375 | S0da=S0ad |
---|
376 | S0db=S0bd |
---|
377 | S0dc=S0cd |
---|
378 | Pdd=2*(Exp(-Xd)-1+Xd)/Xd^2 |
---|
379 | S0dd=Nd*Phid*vd*Pdd |
---|
380 | |
---|
381 | IF(lCASE == 1) |
---|
382 | S0aa=0.000001 |
---|
383 | S0ab=0.000002 |
---|
384 | S0ac=0.000003 |
---|
385 | S0ad=0.000004 |
---|
386 | S0bb=0.000005 |
---|
387 | S0bc=0.000006 |
---|
388 | S0bd=0.000007 |
---|
389 | S0cd=0.000008 |
---|
390 | ENDIF |
---|
391 | |
---|
392 | IF(lCASE == 2) |
---|
393 | S0aa=0.000001 |
---|
394 | S0ab=0.000002 |
---|
395 | S0ac=0.000003 |
---|
396 | S0ad=0.000004 |
---|
397 | S0bb=0.000005 |
---|
398 | S0bc=0.000006 |
---|
399 | S0bd=0.000007 |
---|
400 | ENDIF |
---|
401 | |
---|
402 | IF(lCASE == 3) |
---|
403 | S0aa=0.000001 |
---|
404 | S0ab=0.000002 |
---|
405 | S0ac=0.000003 |
---|
406 | S0ad=0.000004 |
---|
407 | S0bc=0.000005 |
---|
408 | S0bd=0.000006 |
---|
409 | S0cd=0.000007 |
---|
410 | ENDIF |
---|
411 | |
---|
412 | IF(lCASE == 4) |
---|
413 | S0aa=0.000001 |
---|
414 | S0ab=0.000002 |
---|
415 | S0ac=0.000003 |
---|
416 | S0ad=0.000004 |
---|
417 | S0bc=0.000005 |
---|
418 | S0bd=0.000006 |
---|
419 | ENDIF |
---|
420 | |
---|
421 | IF(lCASE == 5) |
---|
422 | S0aa=0.000001 |
---|
423 | S0ab=0.000002 |
---|
424 | S0ac=0.000003 |
---|
425 | S0ad=0.000004 |
---|
426 | ENDIF |
---|
427 | |
---|
428 | IF(lCASE == 6) |
---|
429 | S0ab=0.000001 |
---|
430 | S0ac=0.000002 |
---|
431 | S0ad=0.000003 |
---|
432 | S0bc=0.000004 |
---|
433 | S0bd=0.000005 |
---|
434 | S0cd=0.000006 |
---|
435 | ENDIF |
---|
436 | |
---|
437 | IF(lCASE == 7) |
---|
438 | S0ab=0.000001 |
---|
439 | S0ac=0.000002 |
---|
440 | S0ad=0.000003 |
---|
441 | S0bc=0.000004 |
---|
442 | S0bd=0.000005 |
---|
443 | ENDIF |
---|
444 | |
---|
445 | IF(lCASE == 8) |
---|
446 | S0ab=0.000001 |
---|
447 | S0ac=0.000002 |
---|
448 | S0ad=0.000003 |
---|
449 | ENDIF |
---|
450 | |
---|
451 | IF(lCASE == 9) |
---|
452 | S0ac=0.000001 |
---|
453 | S0ad=0.000002 |
---|
454 | S0bc=0.000003 |
---|
455 | S0bd=0.000004 |
---|
456 | ENDIF |
---|
457 | |
---|
458 | S0ba=S0ab |
---|
459 | S0ca=S0ac |
---|
460 | S0cb=S0bc |
---|
461 | S0da=S0ad |
---|
462 | S0db=S0bd |
---|
463 | S0dc=S0cd |
---|
464 | |
---|
465 | Kaa=0 |
---|
466 | Kbb=0 |
---|
467 | Kcc=0 |
---|
468 | Kdd=0 |
---|
469 | |
---|
470 | Kba=Kab |
---|
471 | Kca=Kac |
---|
472 | Kcb=Kbc |
---|
473 | Kda=Kad |
---|
474 | Kdb=Kbd |
---|
475 | Kdc=Kcd |
---|
476 | |
---|
477 | Zaa=Kaa-Kad-Kad |
---|
478 | Zab=Kab-Kad-Kbd |
---|
479 | Zac=Kac-Kad-Kcd |
---|
480 | Zba=Kba-Kbd-Kad |
---|
481 | Zbb=Kbb-Kbd-Kbd |
---|
482 | Zbc=Kbc-Kbd-Kcd |
---|
483 | Zca=Kca-Kcd-Kad |
---|
484 | Zcb=Kcb-Kcd-Kbd |
---|
485 | Zcc=Kcc-Kcd-Kcd |
---|
486 | |
---|
487 | DenT=(-(S0ac*S0bb*S0ca) + S0ab*S0bc*S0ca + S0ac*S0ba*S0cb - S0aa*S0bc*S0cb - S0ab*S0ba*S0cc + S0aa*S0bb*S0cc) |
---|
488 | |
---|
489 | T11= (-(S0bc*S0cb) + S0bb*S0cc)/DenT |
---|
490 | T12= (S0ac*S0cb - S0ab*S0cc)/DenT |
---|
491 | T13= (-(S0ac*S0bb) + S0ab*S0bc)/DenT |
---|
492 | T21= (S0bc*S0ca - S0ba*S0cc)/DenT |
---|
493 | T22= (-(S0ac*S0ca) + S0aa*S0cc)/DenT |
---|
494 | T23= (S0ac*S0ba - S0aa*S0bc)/DenT |
---|
495 | T31= (-(S0bb*S0ca) + S0ba*S0cb)/DenT |
---|
496 | T32= (S0ab*S0ca - S0aa*S0cb)/DenT |
---|
497 | T33= (-(S0ab*S0ba) + S0aa*S0bb)/DenT |
---|
498 | |
---|
499 | Y1=T11*S0ad+T12*S0bd+T13*S0cd+1 |
---|
500 | Y2=T21*S0ad+T22*S0bd+T23*S0cd+1 |
---|
501 | Y3=T31*S0ad+T32*S0bd+T33*S0cd+1 |
---|
502 | |
---|
503 | X11=Y1*Y1 |
---|
504 | X12=Y1*Y2 |
---|
505 | X13=Y1*Y3 |
---|
506 | X21=Y2*Y1 |
---|
507 | X22=Y2*Y2 |
---|
508 | X23=Y2*Y3 |
---|
509 | X31=Y3*Y1 |
---|
510 | X32=Y3*Y2 |
---|
511 | X33=Y3*Y3 |
---|
512 | |
---|
513 | ZZ=S0ad*(T11*S0ad+T12*S0bd+T13*S0cd)+S0bd*(T21*S0ad+T22*S0bd+T23*S0cd)+S0cd*(T31*S0ad+T32*S0bd+T33*S0cd) |
---|
514 | |
---|
515 | m=1/(S0dd-ZZ) |
---|
516 | |
---|
517 | N11=m*X11+Zaa |
---|
518 | N12=m*X12+Zab |
---|
519 | N13=m*X13+Zac |
---|
520 | N21=m*X21+Zba |
---|
521 | N22=m*X22+Zbb |
---|
522 | N23=m*X23+Zbc |
---|
523 | N31=m*X31+Zca |
---|
524 | N32=m*X32+Zcb |
---|
525 | N33=m*X33+Zcc |
---|
526 | |
---|
527 | M11= N11*S0aa + N12*S0ab + N13*S0ac |
---|
528 | M12= N11*S0ab + N12*S0bb + N13*S0bc |
---|
529 | M13= N11*S0ac + N12*S0bc + N13*S0cc |
---|
530 | M21= N21*S0aa + N22*S0ab + N23*S0ac |
---|
531 | M22= N21*S0ab + N22*S0bb + N23*S0bc |
---|
532 | M23= N21*S0ac + N22*S0bc + N23*S0cc |
---|
533 | M31= N31*S0aa + N32*S0ab + N33*S0ac |
---|
534 | M32= N31*S0ab + N32*S0bb + N33*S0bc |
---|
535 | M33= N31*S0ac + N32*S0bc + N33*S0cc |
---|
536 | |
---|
537 | DenQ1=1+M11-M12*M21+M22+M11*M22-M13*M31-M13*M22*M31 |
---|
538 | DenQ2= M12*M23*M31+M13*M21*M32-M23*M32-M11*M23*M32+M33+M11*M33 |
---|
539 | DenQ3= -M12*M21*M33+M22*M33+M11*M22*M33 |
---|
540 | DenQ=DenQ1+DenQ2+DenQ3 |
---|
541 | |
---|
542 | Q11= (1 + M22-M23*M32 + M33 + M22*M33)/DenQ |
---|
543 | Q12= (-M12 + M13*M32 - M12*M33)/DenQ |
---|
544 | Q13= (-M13 - M13*M22 + M12*M23)/DenQ |
---|
545 | Q21= (-M21 + M23*M31 - M21*M33)/DenQ |
---|
546 | Q22= (1 + M11 - M13*M31 + M33 + M11*M33)/DenQ |
---|
547 | Q23= (M13*M21 - M23 - M11*M23)/DenQ |
---|
548 | Q31= (-M31 - M22*M31 + M21*M32)/DenQ |
---|
549 | Q32= (M12*M31 - M32 - M11*M32)/DenQ |
---|
550 | Q33= (1 + M11 - M12*M21 + M22 + M11*M22)/DenQ |
---|
551 | |
---|
552 | S11= Q11*S0aa + Q21*S0ab + Q31*S0ac |
---|
553 | S12= Q12*S0aa + Q22*S0ab + Q32*S0ac |
---|
554 | S13= Q13*S0aa + Q23*S0ab + Q33*S0ac |
---|
555 | S14=-S11-S12-S13 |
---|
556 | S21= Q11*S0ba + Q21*S0bb + Q31*S0bc |
---|
557 | S22= Q12*S0ba + Q22*S0bb + Q32*S0bc |
---|
558 | S23= Q13*S0ba + Q23*S0bb + Q33*S0bc |
---|
559 | S24=-S21-S22-S23 |
---|
560 | S31= Q11*S0ca + Q21*S0cb + Q31*S0cc |
---|
561 | S32= Q12*S0ca + Q22*S0cb + Q32*S0cc |
---|
562 | S33= Q13*S0ca + Q23*S0cb + Q33*S0cc |
---|
563 | S34=-S31-S32-S33 |
---|
564 | S41=S14 |
---|
565 | S42=S24 |
---|
566 | S43=S34 |
---|
567 | S44=S11+S22+S33+2*S12+2*S13+2*S23 |
---|
568 | |
---|
569 | Nav=6.022045E23 |
---|
570 | Lad=(La/va-Ld/vd)*SQRT(Nav) |
---|
571 | Lbd=(Lb/vb-Ld/vd)*SQRT(Nav) |
---|
572 | Lcd=(Lc/vc-Ld/vd)*SQRT(Nav) |
---|
573 | |
---|
574 | Int=Lad^2*S11+Lbd^2*S22+Lcd^2*S33+2*Lad*Lbd*S12+2*Lbd*Lcd*S23+2*Lad*Lcd*S13 |
---|
575 | Int=Int*scale + background |
---|
576 | |
---|
577 | // print Q,Int,Nd,var[5],lCASE |
---|
578 | |
---|
579 | return Int |
---|
580 | |
---|
581 | End |
---|
582 | |
---|
583 | // this is all there is to the smeared calculation! |
---|
584 | Function SmearedRPAForm(s) :FitFunc |
---|
585 | |
---|
586 | Struct ResSmearAAOStruct &s |
---|
587 | |
---|
588 | ////the name of your unsmeared model is the first argument |
---|
589 | s.yW = Smear_Model_20(RPAForm,s.coefW,s.xW,s.resW) |
---|
590 | |
---|
591 | return(0) |
---|
592 | End |
---|
593 | |
---|
594 | //wrapper to calculate the smeared model as an AAO-Struct |
---|
595 | // fills the struct and calls the ususal function with the STRUCT parameter |
---|
596 | // |
---|
597 | // used only for the dependency, not for fitting |
---|
598 | // |
---|
599 | Function fSmearedRPAForm(coefW,yW,xW) |
---|
600 | Wave coefW,yW,xW |
---|
601 | |
---|
602 | String str = getWavesDataFolder(yW,0) |
---|
603 | String DF="root:"+str+":" |
---|
604 | |
---|
605 | WAVE resW = $(DF+str+"_res") |
---|
606 | |
---|
607 | STRUCT ResSmearAAOStruct fs |
---|
608 | WAVE fs.coefW = coefW |
---|
609 | WAVE fs.yW = yW |
---|
610 | WAVE fs.xW = xW |
---|
611 | WAVE fs.resW = resW |
---|
612 | |
---|
613 | Variable err |
---|
614 | err = SmearedRPAForm(fs) |
---|
615 | |
---|
616 | return (0) |
---|
617 | End |
---|